C SWAN BANDS IN COMETS* 


Nt><. 5” 

RALPH E. STOCKHAUSEN AND DONALD E. OSTERBROCKf • 


The relative populations of the vibrational levels of the X 3 I1 and A 3 n electronic levels are 
calculated assuming the fluorescence mechanism. Pure vibrational transitions are taken into 
account by making a rough estimate of the magnetic-dipole transition probabilities. Roth the 
approximate method, Rosseland’s theory of cycles applied to a three-level molecule, and the 
accurate solution of the equations of statistical equilibrium for a 10-level molecule, give similar 
results. The excitation temperatures derived from these relative populations agree satisfactorily 
with the observations of the Swan bands by McKellar and Climenhaga (1953) for various sun- 
comet distances. Finally, an estimate is made of the infrared radiation, due to pure vibrational 
transitions, expected from a. bright comet. The expected amount of radiation is small and will 
be difficult to detect. 


It is well known that the C 2 Swan bands ob- 
served in comets are excited by the mechanism of 
resonance fluorescence (Stawikowski and Swings 
I960). Observations give excitation temperatures 
(defined by a Boltzmann distribution of the rela- 
tive populations of vibrational levels of the X 3 II 
ground electronic term) in the range T eI ~ 2000- 
3000° at sun-comet distances of 0.48 to 1.40 a.u. 
(McKellar and Climenhaga 1953). Other comet- 
ary molecules such as CN, CH, NH and OH are 
observed to have much lower excitation tempera- 
tures, and this is understood qualitatively to be a 
consequence of spontaneous downward radiative 
transitions (vibrational and pure rotational) with- 
in the ground electronic terms of these molecules 
(Wurm 1936). In fact, quantitative calculations 
using estimated transition probabilities approxi- 
mately match the observed excitation temperature 
(Hunaerts 1953, 1957). However C 2 has no per- 
manent electric dipole moment, and the downward 
radiative transitions within the ground term are 
therefore forbidden, so that large populations 
occur in highly excited vibrational levels (Wurm 
1936). 

Houziaux (1960) has calculated the expected 
populations of successive vibrational levels in C 2 , 
taking account only of resonance fluorescence in 

♦Published as Goddard Space Flight Center document X-614-64~X3t, 

August 1904. 

tWaahbum Observatory, Madison, Wisconsin. 


the Swan bands, and assuming all transition prob- 
abilities within the ground X term to be iden- 
tically zero. However the calculations do not 
agree very well with observation in that the cal- 
culated excitation temperatures are higher than 
the observed excitation temperatures and further- 
more vary considerably from one level to the next. 

In fact, however, magnetic dipole transitions of 
the type A2=+l can occur within the ground 
X 3 n term of C 2 , since all the selection rules for 
this type of transition are fulfilled (Van Vleck 
1934). We can see the approximate effect of 
these transitions by applying Rosseland’s theory 
of cycles to a simplified 3 level molecule, taking 
the lowest levels 1 and 2 as successive vibrational 
levels of X *n, and the high level 3 as any vibra- 
tional level of the upper A 3 n term that is con- 
nected by strong radiative transitions with both 
levels 1 and 2. Then the population ratio is 
(Ambartsumyan 1958). 

Ar 32/ * r +~7^(A3i+A32) 

— — We~ h “2\/kT (1) 

til — h» A 

WAn 33/ * r -f - — 2 '(.4 : S i + A 32) 

•A 32 

=e~ h, 2\/kT„ 

where W is the dilution factor, T is the radiation 
temperature and T IX is the observed temperature 
as defined above. 
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In a crude first approximation Au^An^An 
so if WA 31 exp (—hvn/kT)y>An, that is if the 
transition 2— >1 is strongly forbidden, then T, x 
« T « 5730°. Alternatively, if A 3 i exp ( — h v^/kT) 
»A 2 i, that is 2— »1 can occur, then [IT] exp 
( — hv 2 i/kT)~exp( — hv 2 i/kT tI ). For a comet at 
1 a.u. distance from the sun the geometrical dilu- 
tion factor is JF«5XlO -# , (hv^i/kT) *=0A, and 
this limit gives T tz « 180°. The magnetic dipole 
transition probability is (Condon and Shortley 
1951) 

^(aA/b) 2 , (2) 

where a = 1618 cm -1 is the wave number between 
successive vibrational levels, and (aAfb) is the 
matrix element of the magnetic dipole moment, 
which we can crudely estimate to be 1 Bohr mag- 
neton. We thus roughly estimate Az\(=A tib ) 
= 10 -1 sec — 1 , which is of the same order of m 
nitude as TFA 3 iexp(— hvx/kT), but smaller than 
A 3 iexp( — hvn/kT), since for the strongest transi- 
tions in the Swan bands A 3 i«=5Xl0 8 sec — *. Ac- 
cordingly, we expect the calculated excitation 
temperature to be between the two limits above 
but closer to the first, that is, to be of the order of 
a few thousand degrees. 

Therefore accurate calculations of the statis- 
tical equilibrium of cometary C 2 were made, tak- 
ing into account the lowest 5 vibrational levels of 
the ground X 3 n term and also the lowest 5 vibra- 
tional levels of the upper A 3 II term. 

The relative populations were obtained by 
solving the equations of statistical equilibrium. 
Equations similar to those of Houziaux (1960) 
were used, with the addition of pure vibrational 
transitions within the electronic terno. Only 
those pure vibrational transitions with An = ± 1 
were considered. For completeness, pure vibra- 
tional absorptions of solar infrared radiation were 
included, although they had negligible influence 
on the results. Following are the equations of 
statistical equilibrium, lower electronic term, 
i — 1, 5: 

n<- 1 B { - 1 t pt-i t —n t (B ( £ _i p t ( _i 

10 

-\-B ( £+1 pt (+ i 2 B tk P(*)+n <+ 1 

*-« 

10 

(B< f+i Pi i+i ■(■•dt+i dd - 2n*(Z»* £ p*<+A* £ ) — 0 (3) 
*-« 


upper electronic term, j = 6, 10: 
n j-i j Pi- 1 i~ n )[Bj i- 1 Pi j-i+Aj j- 1 

s 

+Bj y+i pj i+ i-(- 2 (Bj p } +Aj)]+n J+ i 

-i 

(B j+ i jp J+lj +A }+ i ,)-f 2 n B } p j =0 (4) 

=i 

The notation is essentially the same as used by 
Houziaux (1960). 

The same values of A rib were used for all vibra- 
tional transitions of both upper and lower elec- 
tronic terms. For the pure vibrational absorp- 
tions, all upper levels were assumed to absorb at 
X = 5.70m while the lower absorbed at 6.18m (these 
are wavelengths of the t> = 0 to t>=l transitions). 
The solar radiation at these wavelengths was ob- 
tained from Allen (1963). Values of the other 
parameters used in equations (3) and (4) are 
listed in Table 1. The columns headed m and n 
refer to the subscripts used in equations (3) and 
(4). The transition probabilities were computed 
from the overlap integrals and the oscillator 
strength using the same values of these quan- 
tities as Houziaux (1960). The radiation densities 
were derived from the “mean monochromatic in- 
tensities” in Table 1 of Minnaert (1953). The 
transition probabilities and radiation densities 
computed here are slightly different from those of 
Houziaux, but either set of values gives essentially 
the same results. 

The equations of statistical equilibrium were 
solved with the aid of an electronic computer. 
For comparison, we show in Table 2 the results of 
our work for the case with A Hb -Q and a sun-comet 
distance of 0.72 a.u. (which is probably the dis- 
tance used by Houziaux) and Houziaux’s results 
for the same case. We believe that there is some 
error in the latter. Table 3 shows the results for 
various values of A, ib and the sun-comet distance. 
In order to allow for uncertainties in the value of 
calculations were made with the assumed 
values A, ib = 1.0, 0.1 and 0.01 sec — 1 . 

It can be seen that the results of Table 3 have 
the same behavior as the three-level case, i.e., 
decreasing vibrational temperature with increas- 
ing vibrational transition probability. It is also 
seen that transition probabilities between A,< b 
= 0.1 and 1.0 sec — 1 give excitation temperatures 
in the 2000-3000° range observed by McKellar 
and Climenhaga. 
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Table 1. — Parameters Used in xe Equations of Statistical Equilibrium 


V* 

V* 

m 

n 

X 

A ' * 
sec -1 

» ' * 

JL>V V 

erg^cnUsec -2 

p/ /(0.72 a.u.) 
erg cm -3 sec 

0 

0 

6 

1 

5165 

6. 36- 10® 

5. 26* 10 18 

1. 135- 10-‘» 

0 

1 

6 

2 

5635 

1.41* 10® 

1.52* 10 18 

1. 257 

0 

2 

6 

3 

6191 

2. 12* 10 5 

3. 02- 10 17 

1. 433 

1 

0 

7 

1 

4737 

2. 67* 10® 

1. 71* 10 18 

1.020 

1 

1 

7 

2 

5129 

3. 23* 10* 

2. 61* 10 18 

1. 129 

1 

2 

7 

3 

5585 

1.93* 10® 

2. 02* 10 18 

1.263 

1 

3 

7 

4 

6122 

4. 60* 10 s 

6. 34- 10 17 

1.425 

1 

4 

7 

5 

6764 

7. 75- 10* 

1. 44* 10 17 

1.530 

2 

0 

8 

1 i 

4382 

3. 42* 10* 

1.73* 10 17 

0. 741 

2 

1 

8 

2 

4715 

4. 07* 10® 

2. 56- 10 18 

1.009 

2 

2 

8 

3 

5098 

1. 47* 10® 

1. 17* 10 18 

1.096 

2 

3 

8 I 

4 

5540 

1.96* 10® 

2. 02* 10 18 

1. 223 

2 

4 

8 

5 

6060 

6. 68* 10 5 

8. 93* 10 17 

1.418 

3 

1 

9 

2 

4371 

8. 61* 10 5 

4. 32* 10 17 

0. 737 

3 

2 

9 

3 i 

4697 

4. 68* 10® 

2. 92* 10 18 

1.001 

3 

3 

9 

4 

5071 

5. 24* 10 5 

4. 10* 10 17 

1. 119 

3 

4 

9 

5 

5502 

1. 89* 10® 

1. 89* 10 18 

1. 211 

4 

2 

10 

3 

4365 

1. 40* 10® 

6. 98* 10 17 

0. 734 

4 

3 

10 

4 

4685 

4. 98* 10® 

1 

3. 07* 10 18 

0. 995 


THIS PAPER 


Table 2. — Relative Populations 


HOUZIAUX (1960) 


Lower Electronic Term 

Upper Electronic Term 

Lower Electronic Term 

Upper Electronic Term 

V* 

N{v')/N(v* = 0) 

v' 

A r (t>')/AT(t>'=0) 

v w 

N(y”)/N(v n =0) 

v f 

N(v')/N(v" =0) 

1 

0. 69 

0 

9. 4* 10" 8 

1 

0. 95 

0 

1 . o- 10- 7 

2 

0. 48 

1 

6.4 

2 

0. 71 

1 

8. 7* 10~« 

3 

0. 33 

2 

4.3 

3 

0.61 

2 

7.9 

4 

0. 23 

3 

2.9 

4 

0. 37 

3 

4.5 



4 

2.0 



4 

6. 1 


From the relative populations given in Table 3 
it is also possible to estimate the expected infrared 
radiation of cometary C 2 . We will make this cal- 
culation for the (1,0) band of the lower electronic 
term, which occurs at 6.18/i. Both the earth- 
comet and sun-comet distances are taken as 1 a.u. 

According to Wurm (1963), the number of C 2 
molecules in the ground state of a typical bright 
comet is 

N(C) «6X10 S2 molecules 


The luminosity in the (1,0) band is then given by 
the expression 

y N ( c JA.*hv, ib 

where N(v"= 1)/N{v” = 0) is a relative population 
given in Table 3. With the assumption A,n= 
0.1 sec -1 we obtain Z,«9X10 18 ergs/sec, or a flux 
received at the earth, assuming the whole comet 
is in the field of the telescope, of about 3X10"* 
ergs/cm* sec. 
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Table 3 . — Relative Populations and Boltzmann Temperatures 


Lower Level 

2 

Upper Level 

v’ 

N(.v")/N(v'=0) 

V* 

N(v')/N(v'=0) 


R = 

=0.4 a.u. A vi b = 

L0 


1 

0. 39 

25X10* 

1 

0.50 

2 

0. 22 

31 

2 

0. 25 

3 

0. 12 

33 

3 

0. 14 

4 

0. 06 

32 

4 

0.09 



A t?t& =0.1 



1 

0. 65 

53X10* 

1 

0. 62 

2 

0. 44 

56 

2 

0. 42 

3 

0.29 

56 

3 

0. 28 

4 

0.20 

56 

4 

0. 19 



A ,<»=0.01 



1 

0. 69 

63X10* 

1 

0. 68 

2 

0. 48 

62 

2 

0. 45 

3 

0. 33 

61 

3 

0.31 

4 

0. 23 

62 

4 

0. 21 


R 

= 0,7 a.u. A rib = i 

.0 


1 

0.20 

14X10* 

1 

0. 39 

2 

0.09 

19 

2 

0. 13 

3 

0. 04 

21 

3 

0. 06 

4 

0.01 

21 

4 

0. 03 



A =0.1 



1 

0. 56 

40X10* 

1 

0.60 

2 

0. 36 

46 

2 

0. 36 

3 

0. 24 

48 

3 

0. 23 

4 

0. 14 

46 

4 

0. 16 



Arib =0.01 



1 

0. 68 

60X10* 

1 

0. C7 

2 

0.47 

61 

2 

0. 45 

3 

0.32 

60 

3 

0. 30 

4 

; 0.22 

60 

4 

0.20 


R =■ 

1.0 a.u. Arib = 

1.0 


1 

0. 11 

11X10* 

1 

0. 34 

2 

0. 04 

15 

2 

0. 08 

3 

0.01 

16 

3 

0.03 

4 

0.004 

17 

4 

0.01 



A rib =0.1 



1 

0. 47 

31X10* 

1 

0. 55 

2 

0. 28 

37 

2 

0.30 

3 

0. 17 

39 

3 

0.18 

4 

0.09 

38 

4 

0. 12 
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Table 3. — Relative Populations and Boltzmann Temperatures ( Continued ) 


Lower Level 


Upper Level 

V* 

N(v')/N(v'=0) 

v' 

N(v')/N(v'= 0) 



0.01 



1 

0.66 

57X10 2 

1 

0. 66 

2 

0. 45 

58 

2 

0. 43 

3 

0. 30 

58 

3 

0.29 

4 

0. 21 

58 

4 

0.20 



= 1.4 a.u. A vih = 

1.0 


1 

0.06 

8X10* 

1 

0.31 

2 

0. 02 

12 

2 

0. 05 

3 

0.006 

14 

3 

0. 01 

4 

0.002 

14 

4 

0.006 



A*- 0.1 



1 | 

0. 35 

22X10* 

1 

0. 48 

2 

0. 19 

28 

2 

0. 22 

3 

0.11 

31 

3 

0. 12 

4 

0. 05 

30 

4 

0. 08 



A rift *0.01 



1 

0.63 

51X10* 

1 

0. 65 

2 

0. 43 

55 

2 

0.41 

3 

0.29 

55 

3 

0. 28 

4 

0.20 

57 

4 

0. 19 


Rough calculations can ba made as to the possi- 
bility of detecting this band. Assuming that the 
minimum power (i.e. Noise Equivalent Power) 
detectable by an infrared detector is about 1 0 -u 
watts; and that 0.1 of the power incident on the 
telescope reaches the detector, then it would 
require a telescope diameter of several hundred 
inches to observe this band (even above the at- 
mosphere). Clearly, detection of the infrared 
lines of Cj are not to be expected. 

In summary, then, with the estimated pure 
vibrational transition probabilities and the as- 
sumption of the fluorescence mechanism, we have 
been able to derive relative populations of a 10- 
level Cj molecule. These results are in agree- 
ment with those suggested by a simplified 3- 
level molecule and with published observations 
of comets. 

It has also been possible to predict the infrared 
emission expected from a bright comet. This 


emission seems too weak to be observed. Similar 

calculations need to be performed on other comet- 
ary molecules to predict their infrared spectra. 
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EVOLUTION OF O STARS. II. HYDROGEN EXHAUSTION 
GRAVITATIONAL CONTRACTION* 

RICHARD STOTHERS 


The evolution of a star of 30 Me is considered from the end of the stable phase of hydrogen 
burning to the onset of helium burning. Ten models are constructed for the hydrogen-exhaus- 
tion ( E ) phase, and six models for the gravitational contraction ( G ) phase. 

The time scale of the j£-phase is so short (8.8 X10 4 years) that the shell source remains 
peaked at 9 = 0.34 and undergoes little hydrogen depletion. Because radiation pressure remains 
strong in the core, convection does not vanish when the hydrogen content at the center falls to 
zero. Gravitational contraction of the core contributes more to the luminosity than shell burn- 
ing, but since the total luminosity changes little, the structure of the envelope is hardly affected. 
On the H-R diagram, the evolutionary track turns back toward the main sequence when X c =0.03, 
and does not turn away again until the shell source becomes important. 

The G-phase begins when Lh, cow /L <0.001, and lasts 9X10 3 years. Although the tem- 
pi ature in the shell increases, hydrogen depletion remains negligible because of the short time 
scale. However, the shell narrows considerably and its peak moves slightly inward for a while. 
The steep temperature gradient outside the shell c T ses the semiconvective z- ie to move inward 
to g=0.48; the hydrogen discontinuity attains a valu \X =0.1. As the luminosity and radia- 
tion pressure in the shell simultaneously increase, the envelope expands. The shell then behaves 
like a node, since the core continues to contract. The convective region near the center shrinks 
asymptotically to a value q =0.06, but at no time does the core approach an isothermal condition. 
The gravitational energy release is nearly uniform throughout the core in ail phases. At the 
onset of helium burning, Tc^l.SXIO 8 °K and p c = 270 gm/cm 3 . Because of the brightening 
shell source, the stellar radius increases rapidly, bringing the evolutionary track to completion 
of the typical S-shaped curve on the H-R diagram. 

When helium starts to burn, the spectral type is B3, and the star is expected never to return 
to the region of O stars during its active life. 




I. INTRODUCTION 

In massive stars when the central hydrogen con- 
tent falls to a value of 0.03, the whole structure 
undergoes a drastic reorganization. Kushwaha 
(1957) first considered the subsequent early phase 
of hydrogen burning in a shell outside the convec- 
tive core for a star of 10 M©. Reiz (1963) has 
reconsidered and extended this work, but Hayashi 
and Cameron (1962) and Hayashi, Hoshi, and 
Sugimoto (1962) have pursued the study of a star 
of 15.6AT© into the entire hydrogen-exhaustion 
phase and then into the gravitational contraction 
phase preceding helium burning. These phases 


have also been considered for stars of intermediate 
mass, as follows: 3.89M© (Hoyle 1960), 4 Mo 
(Hayashi, Nishida, and Sugimoto 1962; Hayashi, 
Hoshi, and Sugimoto 1962), 5 M© (Polak 1962), 
and 7 M© (Hofmeister, Kippenhahn, and Weigert 
1963). It is the purpose of the present paper to 
follow the evolution of a star of 30Af© during 
these phases to the onset of helium burning at the 
center. The previous hydrogen-burning stage 
has been computed in Paper I (Stothers 1963)-. 

II. ASSUMPTION AND DEFINITIONS 
a) Assumptions 


'Published in The Atirophyticai Journal, 140 , 2 ): 610-523, August The same general assumptions as in Paper I are 
i5. 1064. This psoer and "Evolution ot o stars, i Hydrogen Bum- made, except that nuclear-energy generation also 

mg M were baaed on the author’® doctoral dissertation, Harvard . ... , ,, , 

University, 1963 . occurs in a radiative shell around the core and the 
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time rate of change of the physical variables 
must be taken into account during these fast 
evolutionary phases. 

The opacity of stellar material is again assumed 
to be due only to electron scattering, and the 
equation of state is represented by the sum of the 
perfect gas pressure and radiation pressure. 

b) Notation 

In general, we apply the same notations for 
stellar zones as in Paper I, except here we divide 
the radiative intermediate zone into outer (Zone 
Ilia) and inner (Zone Illb) subzones, separated 
by the radiative, hydrogen-burning shell (denoted 
by a subscript s). 

To designate the successive evolutionary phases 
we introduce H (hydrogen burning), E (hydrogen 
exhaustion), and G (gravitational contraction) 
followed by a number indicating the particular 
model in the phase. 

c) Chemical Composition 

The outer radiative envelope retains the initial 
age-zero composition assumed in Paper I : 

X,=0.70, F,=0.27, £,=0.03, X C no = Z,/2. (1) 

The semiconvective zone constantly adjusts its 
composition to maintain convective neutrality. 
A discontinuity in X marks its interface with the 
inner radiative zone, which retains a constant 
gradient in X, left behind by the retreating, ho- 
mogeneous convective core. The distribution of 
X for the last model of the hydrogen-burning 
phase is given by Figures 3 and 4 of Paper I. In 
the core of that model X e = 0.07. 

The change of hydrogen content as an explicit 
function of time, r, and mass fraction, q, is given 
by 


dX 

dr 

dX 

8 t 


— 4r-("adiative region) 

Ed H 



ndq (convective core), 


0 


( 2 ) 

(3) 


where t H is the energy generation due to hydrogen 
burning and 2?h = 6.0X10 1s erg/gm is the energy 
released per CNO cycle. 


d) Nuclear Energy 

The accurate expression for the rate of energy 
release due to the conversion of hydrogen into 
helium via the full CNO cycle has been given by 
Reeves (1962) in terms of the N u +H l rate: 

c n = 7.94 X 10 27 / u 10U i ( r- — c n oXpT 6 2/3 

\ A CNO/ 

exp ( — 152.3 T 6 ~ U3 ) erg/gm sec (4) 

where 

/u.t= 1 + 1.75 p'l'Tf 21 *, (5) 

0u,i = l+O.OO27 77'* 

-0.0037 TV'*— 0.00007T t , (6) 

Au/A' CN o = 0.99 - 0.00067 T t . (7) 

Here fu.i is the weak electron-screening factor, 
0u,i is the correction term to the zero-energy S 
factor, and Xh/.Ycno has been estimated +om 
tables given by Reeves (1962) for temperatures in 
excess of 2.5X10 7 °K. We are assuming every- 
where the full equilibrium abundance of oxygen, 
which, strictly speaking, is not attained in the 
cooler shell source because of the short time scale 
involved. 


e) Gravitational Energy 

The amount of energy involved during a gravi- 
tational contraction or expansion is expressed as 
the difference between the work done by the pres- 
sure and the change in thermal energy of the gas 
and radiation: 




5+80 1 

3 dr 


In n , 


where 


( 8 ) 


1 


^ rn n t 4 1 P 

7 2iH ,T ~ aT ' V T’ 


1-/3 


(9) 


and the time derivative applies to a given mass 
fraction. As Hayashi and Cameron (1962) have 
noted, the logarithmic arguments in equation (8) 
are independent of position in convective regions 
of the star. 


\ 
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III. BASIC EQUATIONS 

The basic equilibrium equations are given in 
Paper I. For Zones I-IIIa the equations were 
transformed and integrated, and the zones fitted 
to each other in the same way as before. Hence- 
forth we shall call this region simply the envelope. 
It is specified only by the luminosity L (through 
the parameter C ) since, down to the chosen ter- 
mination point q f - 0.500, the composition in the 
radiative region has been fixed by the previous 
models. 

In the core (Zones Illb and IV) we employ q as 
the independent variable and logarithms of r,L(r), 
P, and T as dependent variables. The core solu- 
tions are specified by fi c and T c and are integrated 
from the center out through the shell to q f , beyond 
which any energy sources are negligible. The 
customary series expansions near the center 
(Schwarzschild (1958), p. 114) have been modified 
to include radiation pressure and applied at 
g = 0 002. The step value was taken to be 
Ag=0.001. 

Fitting of the envelope and core is made at q f 
in U, V, and L(r) = L. Since L was expected to 
change little during the present evolutionary 
phases, it was sufficient to integrate a short series 
of envelopes at equally spaced intervals of L and 
then to interpolate values of U f and V f at finer 
intervals of L. 

The change in X(q) as a function of time has 
been given by equations (2) and (3). These may 
be approximated by difference equations, with the 
known run of ( t H /X ) as a function of q from the 
previous model. Letting t refer to quantities of 
the previous model, we write in an obvious 
notation 



(q>q <) (10) 

A In X _ 1 [H / *hY 

At “ " (q<i+ qi )EHj Q IXXJ 

+(f)]rf?=-H e (<fc), ( 9 4>9>0) (11) 

+[S e (94) ~S r (94)]f^S 
<741-94 

(<74 l >9>9«) (12) 


for 94*— 94 small. These approximations are 
slightly more accurate than those used by Hasel- 
grove and Hoyle (1956) and by Hayashi and 
Cameron (1962). Outside the convective core, 
in which X c is known (constant), an iterative pro- 
cedure must be used to calculate X from equation 
(10) or equation (12), since (ta/X) depends im- 
plicitly on X through p. Using X from the pre- 
vious model, we calculate a preliminary value of p 
and hence (tH/X). Insertion of (tH/'X) into 
equations (10) or (12) yields a projected value of 
X. Repetition of this procedure yields a definitive 
value of X correct to the second order. 

The method of producing an evolutionary step 
was varied according to the magnitude of the 
ratio of luminosity supplied by gravitational 
contraction to total luminosity, as follow's. 

a) LJL< 0.1 

The method of Hayashi and Cameron (1962) 
was adopted for this case only, although they seem 
io have used it up to LJL ~ 0.6. Here the part of 
the gravitational energy release inclosed in brack- 
ets in equation (8) is expressed as a quadratic 
function of 9, with coefficients extrapolated from 
the previous models: 

1 IcT 

t> ~2 U 3 ) 

The basic equations, where now t- are 

integrated from the center with input values of 
X et 0 e , T e , a 0 , at, 02, 94*, and the run of X 1 and 
( tn/Xy as a function of 9. At 94, A r is deter- 
mined by equation (11). When a 3elf-coiisistent 
model is obtained, the quantity in brackets in 
equation (8) is computed exactly and improved 
values of a 0 , ai, and o 2 are obtained by least 
squares. Another model is integrated with the 
new values, but this second model is alw'ays 
sufficiently accurate for L e /L< 0.1. 

b) 0.1 <L,/L<0.5 

In this case we specify X e and estimate At by 
extrapolation from the previous modeL. Hence 
*, can be computed exactly from equation (8) with 
the given At. In the course of integrations to find 
the correct $ e and T, it is possible continu My to 
improve At by use of equation (11). 

c) L,/L> 0.5 


A In X 
At 
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As L g / L increases, L becomes more sensitive to 
tg and hence to A r. To cope with the increasing 
difficulty of making Ar converge, we choose a 
parameter that loses importance with increasing 
Lg/L, namely, X c . In this case, we take the time 
step by specifying At explicitly. With a value of 
X e extrapolated from the previous models, it is 
possible to obtain a definitive value for X e in the 
same way as for At in Section III6. Eventually, 
when hydrogen burning ceases entirely in the 
shrinking convective core, we need only specif} A t 4 

IV. HYDROGEN-EXHAUSTION PHASE 

The results of the computations for this phase 
are presented in Table 1. It is curious that the 
stellar radius starts to decrease, causing the star 
to turn back toward the main sequence in the 
H-R diagram, as soon as X c «0.03 for a broad 
range of stellar masses (10<Af/Af©<60), even 
though the mass fractions contained in the con- 
vective core are very different (see Table 2). It 
is clear, however, why the radius will shrink when 


X e becomes small. To maintain the luminosity, 
the central temperature must correspondingly 
increase; but since y. e changes little, R~T e ~ l from 
equation (10) of Paper I. 

For a while hydrogen burning in the core is able, 
alone, to maintain the energy balance (Schwarzs- 
child and Harm 1958). However, eventually an 
accelerated contraction of the core is necessary to 
supplement the deficient Ln )a)rt with a growing 
Lg. When the temperature becomes high enough, 
hydrogen burning in a shell outside the core also 
becomes important. It is interesting to note that 
in the above range of masses the peak of this shell 
source always coincides with the mass fraction of 
the convective core when X c «0.07, i.e., just before 
the radius starts to decrease. Table 2 shows that 
q„ the peak of the shell source, is roughly half of 
go, the mass fraction of the core of the initial main- 
sequence model. Furthermore, g 0 itself increases 
by about 0.15 every time the stellar mass is dou- 
bled ; this is a consequence of the increasing impor- 
tance of radiation pressure (diminishing 0), which 
extends the convective regions outward. 


Table 1. — Evolutionary Models of a Star of SO Mq during Hydrogen-Exhaustion (E) Phase 


Models 



1 

2 

3 

4 

5 

6 

7 

8 

9 

10 

log C 

-3. 272 

-3.266 

-3.264 

-3. 262 


-3. 256 

-3. 252 

-3. 248 

-3. 244 

-3. 239 



+0. 717 

+0.719 




+0. 722 


+0. 725 



<fr 

+0. 550 

+0.547 

+0.546 


+0. 545 

+0. 543 


+0. 531) 


+0. 534 

Xt 

+0.650 

+0.647 

+0. 645 

+0. 644 


+0.640 


+0. 634 


+0. 627 

X, 

+0. 589 

+0.584 


+0. 579 

+0. 578 

+0. 573 

+0. 568 

+0.564 


+0. 554 

<7. 

4-0. 342 

+0. 342 

+0. 342 

+0. 342 

+0. 342 

+0. 342 


+0. 342 


+0. 342 

A 

+0. 614 

+0.608 

+0.605 



+0. 604 




+0. 595 

log T, 

+7. 436 

+7. 458 

+7. 487 

+7. 515 

+7. 542 

+7. 561 

4-7. 571 

4-7. 578 

f 7. 584 

+7. 589 

log ( r./Re ) 

+0. 267 

+0. 240 




+0. 130 



+0. 097 

+0. 089 

Qi 

+0. 328 

40. 320 

+0.316 


+0. 302 

+0. 282 

+0. 258 

+0. 236 

+0.205 

+0. 182 

log X. 

— 1. 523 

— 2. 023 

-2. 523 


-3. 523 

-4. 023 

-4. 523 


-5. 882 

-7.311 

A 

+0.564 

+0. 557 

+0. 554 

+0.555 

+0. 555 

+0. 559 


+0. 567 


+0. 577 

log T , 

+7.665 

+7.696 

+7. 726 


+7. 785 

+7.308 

+7. 823 

+7. 833 

+7. 844 

+7. 854 

log p. 

+0. 701 

+0. 795 

+0.886 



+ 1. 140 

+ 1. 193 

+ 1.229 

+ 1.272 

+ 1.312 

cor#/^- ........ 

+0.999 

+0. 995 




+0.599 


+0. 145 


+0. 002 

Lh, •h«Il/£ 

L,/L 

+0. 001 

+0.001 

+0.004 

+0.004 

+0.011 

+0. 014 
+0.043 


+0.094 
+0. 307 

+0, 141 
+0.540 


+0. 217 
+0.755 

+0. 249 
+0. 749 

log (L/L©) 

+5. 439 

+5. 445 

+5. 447 


+5. 451 

+5. 455 

+5. 459 

+5. 463 

+5. 468 

+5. 472 

log (/?//?©} 

+ 1. 143 

+ 1. 136 

+ 1. 113 



+ 1.057 



+ 1.070 

+ 1.079 

log T, 

+4. 548 

+4.554 

+4.565 

+4. 578 


+4. 595 

+4. 596 


+4. 592 

+4. 589 

r (10 4 years) 

0.00 

+5.91 

+7. 74 

+8. 33 

+8.53 
1 

+P.60 

+8.61 

+8. 67 

+8. 69 

+8. 72 
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It is for this reason, the very distance of the 
hydrogen-rich regions, that the temperature at 
X «0.07 is low and shell burning is consequently 
delayed. (In massive stars, 0 is low not so mueh 
by virtue of high temperature as of low density.) 
In stars of lower mass the point at which X c « 0.07 
lies closer to the center, and shell burning begins 
to contribute significantly to the luminosity before 
the core contraction does. The critical mass for 
this to happen is about '20 M©. If q, is small 
enough that the core mass lies below the Schon- 
berg-Chandrasekhar limit (q « 0. 1 ) , thus becoming 
isothermal or even degenerate, then L , may al- 
ways be small. The approach to isothermality 
may be seen in a plot of log T against q for 5 Mo 
(Polak 1962, Fig. 3), where q, = 0.11. It may be 
inferred for 15.6JW© (Hayashi and Cameron 
1962), where </, = 0.17, since Ln, ,heti>L t until 
X c <0.001. In stars with mass'—'lJW© 
(Schwarzschild and Heiberg 1962), the shell 
source supplies practically all the luminosity out- 
side a very large degenerate core. The critical 
stellar mass for the occurrence of degeneracy in 
the core prior to helium burning is~4Af© 
(Hayashi, Nishida, and Hugimoto 1962). 

At 30A/©,g,»0.1 and the core never even 
approaches an isothermal condition. Thus gravi- 
tational contraction proceeds unchecked to in- 
crease T e in accordance with the virial theorem. 
The subsequent increase in «h reduces X t to less 
than 1 percent. Thereafter the drop in X c out- 
runs the increase in 7V, and c H decreases, ~eing 
supplanted by t,~dT e /dr. However, hydrogen 
burning in the core still supplies over half the 
total luminosity even when X, is as low as 10 -4 . 

The distribution of hydrogen in the cential 
regions may be seen in Figure 1. It is clear that, 



Figure 1. — Depletion of hydrogen as a function of mass 
fraction during the hydrogen-exhaustion (E) and gravi- 
tational-contraction (G) phases. Filled circles represent 
the boundary of the convective core in the E models. 

because of the high central temperature and the 
luminosity requirements, X e begins to decrease 
more and more rapidly with q*. This has the 
effect of necessitating an even sharper increase of 
temperature, which consequently exhausts the 
central hydrogen and ignites a region which is still 
relatively hydrogen-rich but close enough to the 
center to feel the temperature rise. Thus Figure 
1 shows why a shell source develops near <7 = 0.34 
for 303/©. We expect this shell to be initially 
broader in massive stars than in stars of lower 
mass because the mean gradient of hydrogen 
throughout this region is gentler (larger qa—q,', see 


Table 2. — Location of Criticai Interface s in Stars of 'ntermediate and High Mass 
When X.-0.03 


M/MO 


9. 

94 

A’. 

Reference 

5 

0. 21 

0. 11 

0. 10 

0. 74 

Polak (1962) 

10 

. 24 


.07 

.90 

Kushwaha (1957) 

Sakashita, Ono, Hayashi (1959) 

15. G_ 

.42 

. 17 

. 16 

.90 

30 

.00 

0.34 

.33 

.70 

Stothere <1003) 

02.7 

0. 75 


0.41 

0. 75 

Schwawchi'd and Hiirm >58) 
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Table 2). The subsequent growth of the shell 
source is shown in Figure 2. 

In Figure 3 the relative contributions to the 
energy generation may be compared. The mag- 
nitude of tf is very nearly constant throughout the 
core for all stages. 



Figure 2. — Logarithm of the energy generation due to 
hydrogen burring as a function of mass fraction. 
Curves are labeled with the model numbers. 



Figure 3. — Lrgarithm of the energy generation as a func- 
tion of ma&j fraction. Solid curves refer to the contri- 
butions from both hydrogen burning and gravitational 
contraction. The dashed curve refers only to gravita- 
tional contraction. Curves are labeled with the model 
number. 

Shrinkage of the convective core is at first slow 
because remains almost constant. Hence the 
evolution of the center approximates Lane’s law, 
, i~T c *. When the shell source becomes impor- 
tant for A",< 10~V e increases at a slightly greater 
rate than T t % , the only concession to isotherraality 
that the core makes. At this point (Model E7), 
L,/L *0.5 and L n , .h.n/L ‘*‘0.1. It seems that as 
soon as L,/‘L grows much larger, contraction of 
the core begins to be offset by an expansion of the 


envelope (cf. Hayashi and Cameron 1962). This 
expansion is due to an increase of Lu , ,h,n, com- 
pensating for the reduction of /, H , rore . Exactly 
why the envelope expands as L H , ,h e ii increases will 
be made evident in the next section. The result 
on the II — R diagram is that the star takes another 
turn away from the main sequence. 

When the hydrogel, conto rt at the center of the 
star vanishes after Mode! E10, the convective core 
does not disappear. This contrasts with the 
results on stars of lower mass (Hoyle I960; Polak 
1962; Hayashi and Cameron 1962), w .i..ce appar- 
ently it is just the nuclear-energy generation that 
maintains the steep temperature gradient. How- 
ever, in the case of 303/©, 0 is still so low that 
convection persists near the center. 

Tne time scale of this phase is short enough 
(8.8X10 4 years) that the shell source, which re- 
mains peaked at q, = 0.342, undergoes very little 
hydrogen depletion. Because of the struggle of 
the various energy sources to maintain the lumi- 
nosity, L increases but slightly. Thus the struc- 
ture of the envelope remains nearly constant. 

V. GRAVITATIONAL CONTRACTION PH«SE 

This phase is taken to st .ri when L», «,or »/L 
<0.001. The hydrogen content of the convective 
core is negligible, and the sole miolear-energv 
source is the hydrogen-burnin ? shell. In i.'odel 
E9 the rise of the shell source uad terminated the 
upward growth of L,/L, which reached a maxi- 
mum value of 0. V 5. This fraction becomes all the 
more impressive when we compare it with 0.56 for 
a star of 15.63/© (Hayashi and Cameron 1962), 
Table 3 shows that even at the onset of helium 
burning (Model G6) the contribution of gravita- 
tional contraction to the luminosity is still a siz- 
able fraction. This suggests that for stars more 
massive than 303/© hydrogen bum.ng in a shell 
may not become dominant. 

During the shell-burning ph.ise- , the computa- 
tional difficulty of fitting models increases. The 
reason lhs in the shift of the U—V curve to the 
left, as seen in Figure 4. This is caused by the 
developing chemical inhomogeneity, which began 
after the initial main-sequence state, and by the 
shell louree, which produces the leftward loo; on 
the U—V plane. Because of the short time ;cule, 
the amount of hydrogen depletion is ”cry mcll 
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Table 3. — Evolutionary Models of a Star of 30 Mo during Gravitational Contraction (G) Phase 


Models 



1 

2 

3 

4 

5 

6 

log C 

—3. 231 

-3. 216 

—3. 195 

“3. 178 

—3. 165 

—3. 152 



+0. 733 

+0. 739 

+0. 749 

+0. 756 

+0. 762 

+0. 769 

qi~q t 

+ 0.530 

+0. 522 

+0.508 

+0. 498 

+0. 489 

+0. 481 

A* 

+0. 620 

+0.606 

+0. 582 

+0.560 

+0.541 

+0. 522 

X, 

+0. 542 

+0. 523 

+0. 490 

+0. 467 

+0. 445 

+0. 425 



+0. 342 

+0. 342 

+0. 338 

+0. 328 

+0. 328 

+0. 328 

p. 

+0. 589 

+0. 576 

+0.556 

+0.553 

+0.553 

+0. 515 

log T, 

+7.596 

+7. 607 

+7. 626 

+7.654 

+7.668 

+7. 683 

log ( rJRo ) 

+0. 073 

+0. 048 

+0. 001 

-0. 058 

-0. 100 

-0. 142 



+0. 155 

+0. 123 

+0. 093 

+0. 078 

+0.069 

+0. 064 

Pc 

+0. 6b7 

+0.604 

+0.631 

+0. 651 

+0.665 

+0. 674 

log T c 

+7. 873 

+7. 907 

+7. 971 

+8. 035 

+8. 101 

+8. 168 

log Pc 

+ 1.386 

+1. 520 

+1. 760 

+1.992 

+2. 217 

+2. 433 

Lhs *h*n/L 

+0. 307 

+0.391 

+0. 485 

+0. 536 

+0. 561 

+0. 575 

U/h 

+0. 693 

+0.609 

+0. 515 

+0.464 

+0. 439 

+0. 425 

log (L/Lo) 

A 5. 480 

+5. 495 

+5. 516 

+5.533 

+5.546 

+5.560 

log ( R/Ho ) 

+ 1.098 

+1. 144 

+1.239 

+1.353 

+1. 492 

+1.680 

log T. 

+4. 581 

+4. 562 

+4. 520 

+4. 467 

+4. 401 

+4..310 

t (10 4 ye^rs) 

0.00 

+0. 10 

+0.30 

+0.50 

+0.70 

+0.90 


and the shell source remains peaked at roughly 
the same mass fraction. Since the shell provides 
the driving mechanism for the expanding envelope 
while the core is contracting, r at the outer bound- 
ary of the shell remains roughly constant. Thus 



u 


Figure 4. — Evolution of a Btar of 30 Mo in the U-V 
plane during the hydrogen-exhaustion (E) and gravita- 
tional-contraction (G) phases. Curves are labeled with 
the model numbers. Dots and jumps represent fitting 
points. The dashed lines in the solutions represent the 
assumed radiative sone (see Paper I). 


it is here that the change of log r with q is greatest 
(see Fig. 5), that is, U is a minimum, isolating the 
energy-producing regions from the envelope. Un- 
less small increments in q are taken in the numeri- 
cal integrations through the shell, the solutions 
diverge. Since the region behaves like an outer 
boundary (£/,— >0), V tends to blow up for small 
changes of the input parameters, fJ e and T e . 



Figure 5. — Logarithm of the radius as a function of mass 
fraction. Curves are labeled with the model numbers. 
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Hayashi, Hoshi, and Sugimoto (1962) have in- 
terpreted the envelope expansion in terms of the 
U—V plane. For values oi the polytropic index 
3<F> 1, where P = Kp 1+1IN and K and N are con- 
stants, the centrally condensed solutions of a 
polytropic envelope form a leftward loop on the 
U—V plane and converge to a ooint U = 0 and 
V = N + 1 as r— >0. For an electi on-scattering en- 
velope, N is 3 near the surface, decreases inside, 
and then approaches 3 again as r— >0. Hence such 
an envelope approximates the polytropic case 
N= 3, where V— >4 as U —> 0 (see Fig. 4). Interior 
to q„ Z/(r)«Z,H. sheii so „hat the radiative index 
becomes large, since (n + l) r ,d~(l— 0)/*L(r). 
Hence the core approximates roughly the iso- 
thermal case N - a ° , which loops outward in the 
U—V plane much as in Figure 4. 

It remains to see how the shell source actually 
expands the envelope. We have already noted 
that, in order to tap additional nuclear fuel, the 
temperature in the core rises and starts shell burn- 
ing. Consequent to hydrogen exhaustion in the 
core, the temperature in the shell must continue to 
rise so that Lh. sh*u can assume more of the total 
luminosity Durden. As Lh, «h«u increases during 
the evolution, (n+l) r<5<! decreases, steepening the 
temperature gradient beyond q, (see F;g. 6). But 
actually n+ 1 can only decrease until the adiabatic 
value is reached. Therefore, to compensate for 
the increasing luminosity, 0 decreases (see Fig. 7). 
This relative increase in radiation pressure expands 
the envelope. 

To make this point more definite, we look at the 
structure of the envelope, which may be con- 
sidered a polytrope of some characteristic index. 


"w •• 



Figure 6. — Logarithm of the temperature as a function 
of mass fraction. Curves are labeled with the model 
numbers. 



Figure 7. — Ratio of gas pressure to total pressure as a 
function of mass fraction. Curves are labeled with the 
model numbers. 


The effective termination point of the envelope 
must be taken just outside the shell peak, where 
U is a minimum on the U—V plane. Now it may 
be shown generally for a polytrope of any index 
and constant 0=0 e that R^(T e /(} e p c ) m . Elir' - 
inating p, we have I2~( 1 — 0 e ) l/2 / fi c T e . Thus the 
expansion may be considered as fundamentally 
due to the relative increase in radiation pressure, 
as stated above. Further, the same result can be 
obtained by examining equation (8). A decrease 
in 0 causes expansion by making t, negative, and 
an increase causes contraction by making c„ posi- 
tive. The consequences may be seen by compar- 
ing Figures 5-7. During phases of roughly 
constant 0, the change of T will determine the 
change in R. 

If we take the stellar core to be a polytrope, then 
we may similarly explain the changes of shell 
radius, r, (see Tables 1 and 3). By extension, a 
star with m shell sources may be subdivided into 
m+1 zones whose structure is assumed to be 
polytropic. Then, approximately, the changes 
of 0 and T at a given shell should determine the 
changes of r at the next outer shell, unless the 
mass fraction of the shell changes by a drastic 
amount. Qualitative agreement is found between 
this picture of evolution and the accurate pub- 
lished results on the evolution of stars of 4Af© 
and 15.6M© from the initial main sequence 
through carbon burning (Hayashi, Nishida, and 
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Sugimoto 1962; Hayashi and Cameron 1962; 
Hayashi, Hoshi, and Sugimoto 1962). 

In this regard we note that the last model be- 
fore helium burning in a star of 1.3 Mq considered 
by Schwarzschild and Selberg (1962) had 0 4 =0.87. 
Since they neglected radiation pressure and used 
0=1 throughout, it might be expected that the 
radii they derive are too small. Indeed the track 
of their theoretical models on the H — R diagram 
deviates leftward of the observed red-giant branch 
in globular clusters. One reason may be the neg- 
lect of varying 0 in their models; apart from other 
factors, for a closer agreement with observations 
0 should be explicitly included. 

Since the envelope evolves independently of the 
core (each with a fixed mass), the expansion of the 
envelope must lower the local density. Since 0 
is also decreasing in the envelope, the temperature 
near q, actually does not change very much ''-«c 
Fig. 6 ). Beyond the shell, however, it \ 
considerably, because the temperature gradient 
is approaching the adiabatic value. This increase 
of the temperature gradient has two major effects: 
( 1 ) The semiconvective zone moves steadily in- 
ward; by Model G 6 it has reached 7 = 0.48. The 
hydrogen discontinuity grows to X 2 — A's=0.10, 
although fi 2 /M 3 = 0.91 because the actual hydrogen 
content is lower at small q. ( 2 ). The shell source 
becomes narrower, on account of the decreasing 
temperature outside (and the decreasing hydro- 
gen content inside). A curious result, not occur- 
ring in stars of lower mass, is that q, moves slightly 
inward at first because of the rising temperature 
inside and the rather broad shell distribution of 
hydrogen, but thereafter remains stationary as it 
encounters less and less hydrogen. The total 
amount of hydrogen depletion (in the shell) is 
0.0004, and the maximum local change in X is 
0.011 at 7 = 0.329. 

In the core, although 0 e steadily rises, the cen- 
tral convection persists and recedes only to a 
roughly constant 74 = 0.06 at the onset of helium 
burning. Hence convection probably never dis- 
appears completely from the center of very 
massive stars. 

Model G 6 corresponds to the onset of helium 
burning. At this point, Lh*/L** 0.004 and T t is 
150 million degrees. This value of T e is almost 
the same as the corresponding value for a star 
of to. CM© because of the steep temperature 


dependence of the energy-generation formula. i 
However, p c is only 270 gm/cm 3 , lower by a factor 1 

of 10 than in the star of 15.6Af©. But in the j 
heavier star the core mass is four times as large. 

The gravitational contraction phase is very 
short (9X10 3 years), although large changes in 
the stellar radius take place. 

VI. FINAL REMARKS 

Although the computed values of the fractional 
luminosities should not be considered exact (since 
we have neglected in the computations to take ac- 
count of the radiative flux that is absorbed in the 
expansion of the envelope), nevertheless Figure 8 
shows up a brief peak in the luminosity due to 
gravitational contraction, whereas the nuclear- 
energy sources are spread over longer time inter- 
vals. Since the total energy available to the star 
as a result of gravitational contraction is small 
compared with the nuclear-energy store, it is 
expected that the time scale of contraction 
phases will be very short. Indeed, the post- 
main-sequence age at the onset of helium burning 
is 4.9 million years, of which 4.8 million are spent 
burning almost half of the initial hydrogen 
content. 

During this time the interior structure of the 
star undergoes drastic changes, especially in the 



212 2.14 2.16 2.18 2.20 2.22 2.24 2.26 

T (10* Y»«r») 

Figure 8. — The fraction of the total luminosity contributed 
by hydrogen burning in the core (solid curve), gravita- 
tional contraction in the core (dashed curve), and 
hydrogen burning in the shell (dash-dot-curve ) is 
plotted against the time elapsed since the last model of 
the hydrogen-burning phase (Model H4 of Paper I). 
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short (0. 1 X 10 6 years) contraction phase, as shown 
in Figure 9. One complication not occurring in 
stars of low'er mass may arise from the existence 
of the semiconvective zone. That is, the inward- 
moving semiconvection brings hydrogen-enrched 
material toward the shell source, which is simul- 
taneously moving outward. A hydrogen flare-up 
might be expected to occur. However, the results 
of Table 3 show that AX is only about 0.1 at the 
radiative-semiconvective interface. Moreover, 
from the results on a star of 15.6Af© (Hayashi 
and Cameron 1962; Hayashi, Hoshi, and Sugimoto 
1962), the total luminosity is expected to increase 
little during helium burning, so that the structure 
of the envelope (and hence of the semiconvective 
zone) will change very slowly; the outward mo- 
tion of <7, is also rather slow. Thus the star could 
probably accommodate quietly the sharp increase 
in X y despite the narrowness of the shell source. 
How r e/er, it seems unlikely that the inner edge of 
the semiconvective zone will ever penetrate the 
shell before the whole envelope becomes convec- 
tlve during the rightward swing of the evolu- 
tionary track in the H — R diagram toward red* 
carbon-burning models. We need only worry 
about possible hydrogen flaring in stars of mass 
greater than 30A/©. 

The evolution on the theoretical H — R diagram 
is shown in Figure 10 for the phases considered in 



Figure 9. — Evolution of the structural zones from the 
initial main sequence to the onset of helium burning. 



Figure 10. — Evolution of a star of 30 Mo on a plot of 
luminosity versus effective temperature, during the 
phases of hydrogen exhaustion and gravitational con- 
traction. Numbers attached to the models (filled 
circles) represent the age in units of 10 4 years. 

this paper. The ages are given in units of 10 4 
years. At r = 9.07 the spectral class becomes 
B0, and at r = 9.67 it is B3. The star is expected 
never to return to the region of O stars during its 
active life. 
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THE RV TAURI POPULATION* 

RICHARD STOTHERS 


Oil:* understanding of the RV Tauri stars has 
come mainly from photometric and spectrographic 
studies based on several tabulations of these stars 
since 1921. 1-9 Although members of this group 
are too distant ffcr accurate parallax dertermina- 
tion, their absolute magnitudes are known through 
the discovery and identification of some in globu- 
lar clusters. 610 It is the purpose of this paper to 
analyze the spatial distribution of the RV Tauri 
stars in the field and to compare properties of 
these stars with properties of similar stars found 
in globular clusters. 

The data are taken mainly from the General 
Catalogue of Variable Stars 8 ( GCVS ), which lists 
92 stars of this class based on criteria that are 
essentially those of Payne-Gaposchkin, Brenton, 
and Gaposchkind RV Tauri stars are supergiants 
with (1) comparatively stable periodicity of light 
variation, whose amplitude may reach 3 magni- 
tudes, (2) alternating shallow and deep minima, 
which occasionally interchange, (3) double-period 
of 39-150 days, and (4) F to K spectral class, 
which is earliest near maximum light. From the 
GCVS list we disregard those variables designated 
RV? and include only stars with periods in the 
range 39-1 10 days, since this is roughly the range 
of RV Tauri stars in globular clusters. Although 
some types of red variables in the same period 
range have light curves resembling those of RV 
Tauri stars, we assume that those variables with- 
out known spectra are yellow. In some cases 
photographic magnitudes were obtained from other 
sources. 4 •• Since the magnitude of CN Centauri 
is uncertain at minimum, it has been omitted 
from the distance determinations. Preston, 
Krzeminski, Smak, and Williams have shown that 
UU Herculis and V453 Ophiuchi are RV Tauri 

'Published in the Publication* of the Aetronotnical Society of tho 
Pacific. 76 (449): 98-105, April 1964. 


stars. 9 Apart from these two variables, their list 
of RV Tauri stars does not otherwise conflict with 
assignments given by the GCVS, in the period 
range 39-110 days. We have used their magni- 
tudes and period for V 453 Ophiuchi. The total 
number of variables whose distances we may 
determine is 45, of which 25 have known spectra. 

Preston et al 9 have defined three groups of RV 
Tauri stars on spectroscopic criteria: A. — All spec- 
tral features indicate type G or K, except that TiO 
bands may occur during deep light minima. B. — 
The type based on the Ca n lines differs from that 
based on the hydrogen lines; strong CN bands 
occur at and around light minima. C. — The spec- 
trum resembles that group B, except that CN 
bands never appear. Moreover, the stars in 
group A are redder than those in groups B and C 
(which cannot be distinguished photometrically). 
The RV Tauri sts-s may be distinguished from 
the yellow semiregulars on the basis of U,B,V 
colors and strength of hydrogen emission. At the 
end of this paper we shall comment on these 
subdivisions. 

The distribution of RV Tauri stars in longitude 
is far from scraggly, contrary to the results indi- 
cated by the smaller numer of stars used in earlier 
investigations. 7 Indeed, Figure 1 shows a strong 
concentration of stars in the direction of the 
galactic center. Their spatial concentration de- 
pends on a knowledge of distances, which we 
may obtain from globular clusters. 

In general, RV Tauri stars in globular clusters 
are specified by tl ? same criteria as those in the 
field. Small deviations from these criteria have 
been discussed by Rosino 6 and Joy. 8 Other dif- 
ferences will be pointed out below. Using the 
mean photographic magnitudes of the nine certain 
RV Tauri stars in clusters (from the writer’s 
previous paper, 10 hereafter called Paper I), we 
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Figure 1 . — Frequency distribution of the field RV Tauri 
stars in longitude. 

may obtain a period-luminosity relation, if the 
distance moduli of the clusters are known. We 
take the moduli from Hogg’s list, 11 where she 
assumed that = 0.0 for the RR Lyrae stars. 

For M 2, M 5, and M 13, we have used the values 
of M, of the RR Lyrae stars corrected for redden- 
ing listed by Arp 12 and have added C.I. = +0.1. 
For M 22, which closely resembles M 13, 13 - 14 and 
for a Centauri, which has a main sequence turnoff 
at the same magnitude as M 3 15 (also listed by 
Arp), we assume analogous values of M pt . For 
M 56 and NGC 6712 we adopt the mean of all the 
known values, M p0 - +0.3. Applying these cor- 
rections to Hogg’s distance moduli, we obtain 
absolute magnitudes of the cluster RV Tauri stars. 
The period-luminosity relation is then <M VQ > = 
— 4.0+0.026P + 0.36, where P is in days, in rea- 
sonable agreement with previous determinations 
which depended on indirect or statistical methods. 

4,6.7, 17. 18, ID, 20, SI 

To calculate the distance of the RV Tauri stars 
in the field, we use 

<M pt > = <w H > — 10—5 log r-a 

where r is in kpc and o = 0.7r if |6 n | <20° or 
0.25 esc | b 11 1 if 1 5 11 1 > 20°. The first formula for 
extinction is used by Payne-Gaposchkin 13 and 
also happens to be the mean of extinctions derived 
for five RV Tauri stars by Kameny. T 

Figure 2 shows the RV Tauri stars projected 
on the galactic plane. Variables with known 
spectra are indicated by open circles, those with 
unknown spectra by filled circles; the crosses 
denote the positions of the sun and the galactic 



Figure 2. — Field RV Tauri stars with known spectra 
(open circles) and unknown spectra (filled circles) are 
projected on the galactic plane. The positions of the 
sun and galactic center are denoted by + signs. The 
circle around the sun has a radius of 1 kpc. 

center. It is difficult to corroborate the previ- 
ously suspected 7 groupings near Ophiuchus (30°), 
Aquila (40°), and Gemini (180°), on the meager 
data available and the uncertain extinction. 
Other dumpings at 60° and 75° are also suspect, 
although the Sagittarius group (0°; is probably 
real. There is no apparent correlation of period 
with position in the Galaxy. 

Figure 3 shows a vertical cross-section of the 
Galaxy, where the sun is assumed to lie at 10 kpc 
from the galactic center. As indicated also by 
Figures 1 and 2, the number of variables increases 
in the direction of the galactic center. Their 
apparent absence at low galactic heights is due to 
interstellar extinction. The distribution strongly 
resembles that of the long-period variables with 
P<250 days, that is, intermediate between the 
distribution of long-period variables with P> 250 
days (Population I) and the distribution of the 
halo RR Lyrae stars (Population II). 13 The 
vertical density gradient y, where log N = x — y 
(r sin b) and N is the relativ e numb er of stars, 
is 1.1, and the median height r sin b is 0.45 kpc. 
Both these quantities are similar to those for the 
long-period variables with Pv, 250 days. 

Perepelkina suspected that the RV Tauri stars 
form an intermediate subsystem. 17 Their rela- 
tively low radial velocities and latitudes indicate 
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Figure 3. — Vertical cross-section of the Galaxy, showing 

the distribution of field RV Tauri stars on the assump- 
tion that the sun lies at 10 kpc from the galactic center. 

Same symbols as in Figure 2. 

a weak Population I, but their excessive blueness 
(like the W Virginis stars) and discontinuous 
velocity curves indicate Population II. 7 The 
yellow semiregulars, however, seem to be real 
Population II objects. 

The general similarity of the RV Tauri stars in 
globular clusters to those in the field and the 
approximately equal ratios of number of RV 
Tauri stars to number of RR Lyrae stars 7 indi- 
cate that field and cluster RV Tauri stars belong 
to the same family. At maximum light the RV 
Tauri stars are usually the brightest cluster mem- 
bers and are found close to the cluster center 
(compared with the long-period variables). 10 
Therefore, because of the long relaxation times of 
globular clusters, it seems unlikely that the RV 
Tauri stars in the field have evaporated from 
clusters, even apart from the observation that 
they do not form a halo population. 

Besides their spatial distribution, some other 
differences between the field and cluster variables 
should be mentioned. While Rosine found no 
period-spectrum relation for the variables in the 
field, 5 Paper I showed such a relation for the 
cluster variables, although the stars seem to be 
separated chiefly into a GO- and a 90-day group. 
Kameny 7 found that the short-period RV Tauri 
stars in the clusters were slightly bluer than those 
in the field. No very long-period RV Tauri stars 


have yet been identified in globular clusters. 
Moreover, in the clusters there is a marked ab- 
sence of RV Tauri stari stars in the period range 
68-87 days. We note in Table I an actual anti- 
correlation of period frequency with the corre- yj 
sponding field vai iables. Lastly, we mention that 


Table I . — Period Frequencies of the RV Tauri Stars 


p 

(clays) 

Globular Clusters 

Field 

RV 

All Variables 

RV (with Sp.) 

RV (all) 

28-47... J 

0 

1 

3 

5 

48-67.... 1 

5 

7 

6 

15 

68-87... J 

0 

2 

15 

23 

88-107... 

4 

13 

3 

5 

>107.... 

0 

16 

5 

24 


no RV Tauri stars in globular clusters have light 
amplitudes exceeding 2 magnitudes, whereas 
about one-third of the field stars in the same 
period range have greater amplitudes. All these 
differences may be due to different ages and 
chemical compositions. 

The results of work by Preston et al .* indicate 
that, at least spectroscopically, the RV Tauri 
stars do not comprise a homogeneous class. Their 
group A of these stars shows TiO at minimum and 
resembles kinematically an intermediate (disk) 
population.* Group B is probably also related 
kinematically to the disk, and we note that it, 
too, shows the presence of metals in the strong 
CN bands. The radial velocities of group C stars 
are very high and suggest a halo, population; the 
absence of CN ia also suggestive. We note that 
the two variables from globular clusters that 
Preston et al. studied probably belong to the 
“halo” group C; it is interesting that, whereas 
groups A and B show a wide range in period, the 
three field variables from group C have periods 
lying in the “forbidden” range of the clusters. 

If the RV Tauri stars plotted on the vertical 
cross-section of the Galaxy (Fig. 3) are distin- 
guished according to groups, then it appears that 
groups B and C adhere very well to the above 
galactic assignments, whereas group A populates 
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both the disk and the halo. However, it is clear 
that RV Tauri stars do not really range far enough 
to separate groups into halo and disk populations. 
Therefore the suggestion by Preston et al. that 
groups B and C belong to a larger family raises 
the further speculation that all the groups form 
really one family. At least this has been indicated 
by previous work. 

In this connection, we note that there is no 
apparent relation between spectroscopic criteria 
and shape of light curve or period. Moreover, 
the colors appear to form a continuous sequence, 
despite the spectroscopic differences. Although 
Preston c' al. could not obtain reliable luminosity 
classes at the dispersion they used, Rosino’s data 5 
suggest that the luminosity class falls with in- 
creasing period, in agreement with the results on 
variables in globular clusters. 10 Comparison of 
the general period-luminosity relations for the 
variables in globular clusters and those in the 
field suggests that any difference in the luminosity 
between the cluster and field RV Tauri stars does 
not exceed 1 mag. (Paper I, Fig. It). We believe 
that the statistical identification of luminosities 
in the clusters and in the field is sufficient to 
outline the general spatial distribution of RV 
Tauri stars, as cot: idered in this paper. 

In conclusion, th° RV Tauri stars seems to 
form an intermediate subsystem between Popu- 
lation I and II. Except for the underabundance 
of metals and the period-frequency anomaly, the 
RV Tauri stare in globular dusters resemble their 
counterparts in the field. Perhaps long-period 
variables in clusters with late integrated spectra 
should be investigated for RV Ta in behavior: in 


this connection we note the 105-day variable in 
NGC 0712 (GA). Certainly the known semiregular 
and irregular variables with undetermined periods 
should be checked for periods in the range 08-87 
days. 
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LONG-PERIOD CONTRIBUTIONS TO THE DISTURB- 
ING FUNCTIONS OF THE EARTH FROM THE 
SEVENTH, NINTH, AND ELEVENTH ZONAL HAR- 
MONICS* 


T. L. l'ELSENTREGEIt AND W. J. WICKLESS 

INTRODUCTION 

It is the purpose of this paper to present explicit 
formulas for the long-period terms due to the 
seventh, ninth, and eleventh zonal harmonics in 
the disturbing function of the earth in the case of 
an artificial earth satellite. The formulas are 
given for terms of the satellite’s orbital elements 
and the Delaunay variables. G. Giaeaglia (1) 
has given general expressions for the long-period 
terms due to any of the zonal harmonics, which 
can be expressed in terms of the orbital elements. 

The apparent differences between the results of 

where 

H = GM 

G = gravitational constant 
Af = mass of the earth 
R = radius of the earth 
J n *= zonal harmonic coefficients (n+2, 3,...) 

P,- Legendre polynomials (n~2, 3,...) 

* = geocentric latitude. 

Here, the earth’s radius is adopted as the unit of length. The seventh, ninth, and eleventh Legendre 
polynomials are 

Pi (sin *) »~(429 sin 7 * -693 sin**+315 sin 1 # -35 sin *) 

ID 

P t (sin «) -—(12155 sin»* -25740 sin V+ 18018 sin** -4620 sin**+315 sin *) 

P„ (sin *) - — ^(88179 sin"* -230945 sinV+218790 sin 7 * -90090 sin»*+ 15015 oin»* -693 sin *). 

'Published M Goddard Spat* fliphl Ctnltr Documnl X-UT-9i-tSS, Augllat 19S4. 
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this paper and those in Giacaglia’s have been veri- 
fied i-.s due to the errors in the latter as it appears 
in the A.J. 

' The contributions of the long-period terms to 
the mean motion of the argument of perigee are 
also given. 

THE DISTURBING FUNCTION 

The eaith’s gravitational potential at a distance 
r from the center of the earth is 

t/*;|i (6in 
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Let 

a - semi-major axis of satellite’s orbit 
e~ eccentricity of orbit 

i = inclination of orbital plane to equatorial plane 
t = mean anomaly 
/=true anomaly 
g = argument o i perigee. 

The Delaunay variables L, G, and H are 

L= y/fia 
G = Ly/l—e* 

H~-G cos i . 

Use is also made of the relations 

sin <f>= sin i sin (f+g) 

r == §*( 1+e cos fi‘ 


The long-period terms in the expansion of U as a Fourier series in l and g are given by 


making use of the relation 



(see reference 2) 




Denoting the long-period parts of Ui,U», and f/n by A 7 F 2 p,A # F 2 p, and AnF 2p , respectively, we have 


Aj F 2 P = 


A9F 2 p = 


21 n*Jie sin i 
16384L 


3*i "J*e sin i 
524288 LW 


[ 10^5 -135^+495^ -429^J ^33 -3og+5^sin g 
-15(3 -69^+209^ -143^) (ll -Hg+sg^sin 3 g 
+33(l -15|+27g -13g) (l -2g+g)sin 5*] 
[210^7 -308^+2002^ -4004^+243 lg)x 


X^7i.. -100lg+385g -35g^sin g -10780^1 -40^+234^ -416— +22l|-^X 

X^39 -65g+29g -3g)sin 30+12012^1 -32^+146^ -200^+85^ X 

/ G 2 G A G*\ / H 2 H A //• //A 

X( 3 -H^+7^ -H-.jsin 5, —2145(l -52g;+ 17|-,)x 


x(> -f,+$ ■ 
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1 1 i.i3 T„p «in iT / H* H 4 H * 7/8 //io\ 

***- ~6 m88^L ^[ 126 ( 21 - 1365 ^+ 13650 ^ -46410^+29393^ -23933^jx 

/ G 2 G* vf G 8 \ 

xf 41423 -77292^+45738^ -8652^+63^-Jsin g 

( H 2 m h 6 //* // l0 \ / G 2 G 4 G® G 8 \ 

1 -61^+570^5 -1802g5+2261^- 8 -969^1 X( 323 -68(^+462-^ -112^+7-^Jsin 3g 

+7020^5 —265^+2130^5 -5746^+6137^ -2261^X^323 -816^+678^ -200^-t 15^sin 5? 

( f/ 2 77 4 77® 77 8 H 10 \ / G 2 G 4 G® G 8 \ 

1 -4lg5+250gj -514^+437^ -133~jx(u -60^+66^ -28^ -3^Jsin Ig 

( H 2 H 8 77® 77 8 77 10 \ / G 2 G 4 G® G 8 \ 1 

1 -25^+90^ -130^+85^ - 21-^5 jx^l “4^+6^ ~ 4 ^®"*'I 5 ) sin 9fll J 


CONTRIBUTIONS TO dg/dt 


Since tha Delaunay set of variables is canonical with respect to the Hamiltonian F, which includes 
A 7 F 2 „ A 2 F 2p , and An F? p , we have 


dg 

dt 


dF 

dG 


(see reference 2). 


Therefore, a computation of d(AiF 2p )/dg(i = 7,9, 11) provides the long-period terms in dg/dt due to the 
seventh, ninth, and eleventh zonal harmonics. 

The results are 


d(A 2 F 2p) _ 21 ft^J' 

~~dG 16384L 3 G 14 c 


H 2 77 4 IP' 
- 135 ^+ 49 % -+ 29^5 


^)(33-3(^+5g) 


r-T— :{l0| (sin 2 i -e 2 )(b -1 
+e 2 sin 2 ^65 -2025^+8415^ -8151^ ^33 -30^+5^ 

+20e 2 sin 2 ^5 -135^+495^ -429^^3 Jsin g 

r / IP H 4 /7®\ / G 2 G 4 \ / H 2 // 4 

-15 (sin 2 i -e 2 )h -69^5+209^7 -143^ Mil -14^+3^ J+e 2 sin 2 rf 39 -1035^5+3553^5 

-2717^ (ll -l/-~+3~j+4e 2 sin 2 ^3 -69^+209^ -143^^7 -3^ jsin 3 g 

+33 [ (sin 2 i -e 2 )(l -15^ +27^ -13^ (l -2^+^+c 2 sin 2 ?|l3 -225^+459^ 

-H) 0 - 2 r>g)^ - <> *} 
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(7« GSX / //2 J/4 //6 /78 //10\ 

-112^+7^)++ sin 2 ; / 21 -1403^+14250^ -48654^+65569^ -30039^1 

/ G 2 G 4 G® G 8 \ ( II 2 II 4 //« //* 

Xl 323 -680^+462^ -112— +7^,) +8+ sin 2 ill -61^+570^ -1302^+2261— 

77 1# \ G 2 / G 2 G 4 G 6 \1 

-969^jx r ^170 -231^+84^ -7^ jjsin 3«7 

+7020^ (sin 2 i-c 2 )(?> - 205^ 2 + 2130^- 5740^+6137-^,- 2261 ^X^323 - 816^+678^ 

G 6 G 8 \ / 772 77 4 77 6 77 8 77 10 \ 

-200^^+15^J+e 2 sir. 2 ?( 105 -6095^+53250^ -155142^+177973^ -70091 

( G 2 G 4 G 6 G 8 \ / 77 2 774 77« 77 s 

323 -816^+678^ -200^+15^+24+ sin 2 ii 5 -265^+2130^ -5746^+6137^ 

7 /io\ ( 72 / (72 (74 

-2261^X^68 -113^+50^ -5^)Jsin 5 g 

[ ( H 2 77 4 77« 77 8 77 10 \ / G 2 G 4 G® G 8 \ 

(sin 2 i —e 2 )ll -41^+250^ -514^+437^ -133^jx(l9 -60^+66^ -28^ 6 +3^j 

( 77 2 77 4 77 6 77 8 77 10 \ / G 2 G 4 G® G 8 \ 

21 -943^+6250^ -13878^+12673^ -4123^ i5 jxf 19 -60^+66^ -28^+3^ J 

+24e sin 2 ?^1 — 41^+250^j —514^+437^ -133^yX^^5 — 11 L 2 ^^L 4 ~L 8 / j Sm ' 

[ ( H 2 77 4 77 6 77 8 77 10 \ / G 2 G 4 G 6 G 8 \ 

(sin 2 1 -+)(^1 -25^+90^ -130^+85^ -2% 0 jx(^l 

( H 2 H* 7/6 77» //io\ / (72 (74 (7# (js \ 

21 -575^+2250^ -3510^+2465^ -651^jx^l -4^+6^ -4^+^ J 

( 7/ 2 77 4 77® 77 8 77 10 \ G 2 / G 2 G 4 G®\1 ) 

1 -25^+90^ -130^+85^ -21-^Jx^l -3^+3^ -^jjsin 9 gj. 
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NOTES ON VON ZEIPEL’S METHOD* 


GIORGIO E. O. GIACAGLIA 


1. INTRODUCTION 

Since the rediscovery of von Zeipel’s method by 
D. Brouwer (1959) and its successful application 
to the problem of artificial satellites, many other 
problems have been solved by that same method, 
'.'.ms proving its great applicability. It is the pur- 
pose of these notes to present the general equa- 
tions of von ZeipePs method and discuss briefly 
their applicability. 

The reduction of the order of a differential 
canonical system can, in theory, be performed by 
obtaining, one way or another, integrals of the 
system. One of them is the Hamiltonian itself 
when it is time independent. Actually, this in- 
tegral of the system (physically its “energy”), can 
describe completely the geometry of the solutions 
in a phase space of 2 n dimensions where 2 n is the 
order of the system. When this order is 2, then 
the solution is completely specified and the use of 
the Hamiltonian reduces it to a first order differ- 
ential equation which can be integrated by quad- 
rature. The introduction of p integrals in a sys- 
tem of n degrees of freedom (2n w order), reduces 
it to one otn—p degrees of freedom which can be 
integrated immediately when n— p< 1 (where, 
of course, p cannot be greater than n). 

A few comments can be made with respect to 
the more famous methods of reduction to show 
their eventual connection with von ZeipePs 
method. 

2. FROM HAMILTON TO VON ZEIPEL 

In the discussion that follows only methods 
that have been used in connection with differen- 
tial systems describing the motion of a physical 
system are considered. The presentation does 
not necessarily follow a chronological order. 

•Published aa Goddard Space Flight Center Document X-6t7-43A-t6t % 
June 1964. 


Consider then a system of n degrees of freedom 
given by 2 n first order canonical equation 

. dH 

X> dy, 

(j = l, 2, •••, n) (1) 

dH 

V * dx, 

where the Hamiltonian H = (xi,...,:r n ,f/ 2 ,... 2 /n) 

is presumed to be time independent. If this is 
not the case, the introduction of time as a new 
canonical coordinate x n +i (the associated momen- 
tum being — H ) always reduces the latter to the 
former case. The degree of freedom will however 
increase by one. 

A canonical transformation of ihe variables 
(x,y) to new variables (x',y f ) will be, in this ex- 
position, equivalent to the problem of finding a 
generating function S-S{x f ,y,t) such that 


ax j 

(j = 1, 2, ») (2) 

dS 

x * = w; 

It is easily seen that this is a sufficient condition 
to satisfy the Jacobi-Poincare relation 

iixjdyj-x'jdy'^dW. (3) 

which is valid whether or not S is an explicit func- 
tion of time. The Hamiltonian of the new sys- 
tem will be equal to that of the old one inasmuch 
as one is obtained from the other by introducing 
the transformation of variables expressed by 
Equations (2) when dS/dt=0. 

a. HA MIL TON-JA COB I— The method intro- 
duced by Hamilton and Jacobi consists in obtain- 
ing a canonical transformation such that the new 
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Hamiltonian is identically zero. In such a case, 
the new variables are all constants. 

b. LINDSTEDT’S METHOD— Lindstedt’s 

method is a particular application of the Hamil- 
ton-Jacobi method for cases where the Hamilto- 
nian is expanded in terms of small parameters. 
In this particular case the solution gives the coor- 
dinates as linear functions of time and the mo- 
menta as constants (usually called action angle 
variables). The comparison with the Hamilton- 
Jacobi method is purely formal since the method 
devised by Lindstedt- is quite original. Actually, 
the real difference between von Zeipel’s and this 
method is that Lindstedt does not make use of a 
generating function. 

c. WHITTAKER’S METHOD (solution by 

series). This method obtains n integrals of the 
system by reducing the Hamiltonian to a function 
of the products In this 

case, since 


, dH dH 

Xj dyj dp Xj 


Vr 


dxj 


dH 


dpj 


Vi, 


it follows that 

ijili+X)y]=0 or p, = const 0'+ l,2,...,n). 

d. DELAUNAY’S METHOD — This method, 
as Lindstedt’s, can be applied only when the 
Hamiltonian consists of a “zero order” part (the 
corresponding system having a known solution) 
and a “disturbing function” that has a small nu- 
merical factor. The basic approach of the von 
Zeipel’s method is the same as that of Delaunay’s 
method; however, the latter one makes no use of a 
generating function and breaks the disturbing 
function into parts which are treated separately. 
The Hamiltonian must be constructed after the 
transformation is performed for each particular 
part. 

A few more techniques could be mentioned but 
one deserves more attention than all the others. 
The concept of adiabatic invariants in Quantum 
Mechanics is quite analogous to the concept of 
“mean variables” in von Zeipel’s method, or to a 
certain extent to what Whittaker calls Adelphic 
Integrals. 


3. THE VON ZEIPEL'S METHOD (1916) 

It has been quite common, after Delaunay, to 
use the negative of the Hamiltonian. Thus, if 
F- —H and if and L,(j= 1,2, ...,n) 

are the coordinates and momenta respectively, 
then 


0=1, 2, •••, n) (4) 



suppose 

F = F(l,L-,*) (5) 

where e is a “small parameter” and l and L indi- 
cate the sets {( h ...,(n) and A 

canonical transformation involving the parameter 
e will be given by a generating function 

S=S((,L*; e) 

such that 


T dS 

(j= 1. 2, «) 

__ dS 

J ~dL* ( 6 ) 

where ( t*,L *) are the new coordinates and mo- 
menta. If F* is the negative of the new Hamil- 
tonian, then we assume 

F*((*,L*; «) «); «) (7) 

or, from Equations (6), 

In a more restrictive sense it is assumed that the 
series 

F= JtfFM, L ) (9) 

*- 0 

represents the negative of the Hamiltonian to the 
required degree of precision and converges to 
F{(,L;t) as N— > <» . From this point F is written 
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as F without danger of confusion. Furthermore, are convergent for sufficiently small e. 


it is assumed that 


S(t, L*; t) 
F k ((, L) 

L*;t ) 


By the conservation property 

*«*4l 


+4 

k - 0 


Jfc-0 


V de* A-o (11) 


where it is important to note that dS/dl contains 
€ through Equation^ (6). Equating the coeffi- 
are developable in Taylor’s series in the neighbor- cients of like powers in t in both sides of Equation 

(11) gives a system of partial differential equa- 
tions in S and F*. The next step is obtaining 
this system. 


hood of «=0, so that the series 
5= 2 t k S k ((, L*)_ 0 

<io > 


*- 0 


'««0 


4. INFERENTIAL EQUATIONS OF THE VON 
ZEIPEL’S METHOD 

The m ,h derivative of F k with respect to e at the 
point e = 0 is obtained as follows. 


Consider 


dF\ 

dt 


t A dF k dLi _ ^ dF \ d / d>S\ 

2->dLi dt ^jdLi ^k\dti) 

Using Equation (10) it follow j that 

(*&-**£ a(v £Ml 
\ dt ) 4dLj dt 4* dli 


L *) v 

= 4 


dL 


<- 1 i - 1 


*/ dS A 

i\d({ /«_o 


Let us now compute 


where 


dt m ~\ dLj 


dt 

F k 


Applying Leibniz’ formula, this becomes 


d m ~ l ( * nin(m Ai' ''“/ni-hdV' 1 d m - l -’(dF k \ 

dt m -\ dl ( ) 4 \ v ) dt’ dt m - l -\dL t )' 


For < = 0 the only possible choice is j<m . Then 


if, -r _1 w-i)'f— (-^)1 

«-o \ «-0 


d m ~ l i 


dt m ~ »* 



Now using Equation (12) 


m -i (f) 

«-0 <-l 1 \ / «-0 •— 0 


It is now desirable to rewrite the above equation as 


(d n F k \ V V (m-l)\( h \ d m ~ n (d l F k \ 

\ ^ dfJ A-'UlJ • 

«“0 <j«i ' f <»oL > / J«~o 


( 12 ) 


( 13 ) 
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Equation (13) is now applied to find 


a- 

3 

1 

(M*\ 


\dLj 


<-o 


The result is 


s. 

if d 2 F„ \ 

< 

1 

<x> 

dt"-> {-’ t \dL h dLj 


d* n -’\dL t J 


The process is repeated up to the point where 

m—ji —jt — • • • — j n = 0, 

so that 

' <*— 'r'» } n_( aX \ 1 _/ e s F k \ 

d« m- VV’ * ’ '~ t N\dLt l dLi i ' • *dL< / \dL<.* • *dZ/< / 

Substituting these successive derivatives into Equation (13), it follows that 
(d m F K \ A A. (m-l)l (dSn\ A (m-j.-l)! ( dSj\ 

V o't m ^ ‘^r 1 (m-j 1 )!\dA 1 /«-o *4 (m-ji-jiV- \ d( h J,-o 


i—0 


(14) 


<i-> 'i- 1 


<2—1 1 


( m ~ji -jt ~j ») !\ <3A s /«-o ^ ^ (m-ji-ji /at)! \ 5Aat/ 

x( \ 

\dLi x dLi 2 m • •3L<^/«-o 




The numerical factors are reduced to 


rn! 


and the second summation does not run in general up to infinity but to a limit given by condition (14). 
Thus the above relation becomes 


(Sp'X.-s n (i)-tc(wA n (^) Xv- f kX 

(m) p-1 \t>-l/ p-1 W 


(15) 


where 2 stands for summation over all possible positive integers j, whose sum is m (according to Equa- 
ls) N n 

tion (14)), and the first product n refers to the summation signs 2. There are n of these integers. 

p-i <»- i 

Equation (15) will be valid even for m*=0 with the definitions 


d°F* 


■F* 


and 


dt° 

C(0;-)«1. 
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From Equations (10) and (11) it follows that 

X* * ’ X* Xm!\d« m /,_ 0 X X m! \ d« m /,_o X Xw!\ dt m /,»o ‘ 

— » ' ■— » n r-0 m-0 ' 


*— 0 *— 0 m— 0 ’ ' *— 0 77t—0 

The substitution of these results into Equation (14) leads to 


'-Sisn ($)«*•« *• -'«n /• 


I (m) p— 1 

In a complete similar way, if 


*- 0 ' 


then 


F* 


-22 xn(.iW*. 

r-o m-0 (m) p-l Vp-l/ p-1 N ' 4-0 > i x itf' 4-0 


(16) 


(17) 


It is important to note that in Equation (16), «=0 is equivalent to L r = dS 0 / dt r (r = 1,2, • • *,n), and in 
Equation (17) «=0 is equivalent to C* = dSa/dL*(r= 1,2, • •*,«), according to Equation (10). The 
equality of factors of the same power of c in Equations (16) and (17) gives the partial differential equa- 
tions for the von Zeipel’s method 


X XnfxW;*, 



Sd h 


dL* 



for i> = 0, 1, 2, • • • . 

For instance, Equation (If) gives: 
* = 0 



(18) 


(19) 


( 20 ) 


JdSo , A , A/ as, dFV\ , A/ as, *£?\ .1 V / dSi d'Ft \ 

*\dL* ,L> J^Z\dL? dl* }„ dS 0 ^±\dL! dl* ) dSo 2 4 \dL* dL * dlTdl*J„ 


dS<t 

w* 


( 21 ) 
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^ + *>(** **i\ x. ^ 1/dSi dSi W\ \ 

\ dL t ) . _3<f? 0 4\d£ ( dL t ) T dSo 4 2\d^, dLidLj) r dS 0 

' L ‘-% -‘ 1,_ »S l ‘~sT," 

+ V i/^1 AS. _»t_\ -A /M AS. w, \ , A/ax. af-A 

^ 6\ d( t dtj dt k dLidLfdLk) j _ d<So H- ^ \ <3v ( d(j dL t dLj)j _ dS 0 4\df< dLif . dS Q 


r j v*_7|)-r — > i 


r T _ \JK1Q 1 

L * - ^* ‘- 1 


— F*( J^*2 T * \ 1 y?( d*'* 1 xv/ ^ dF 2 \ ^ l/ <?^ i <W>i d 2 /^* \ 

\dL*’ rfAdL* dt*) dSo + MdL* dft) dS 0 + 4 2 \dLf dZ* <LS 0 

,_1 *‘~dZ? ‘- 1 {i ~lL* utml • ^ = 1Z* 

, l/ dfli diS't \ .4 / d£i 3^2 d;Fo \ 

* 6\3L* dL* dL* d^ -1- 4 VLfaZ* d(*d(f) dSp 


r j* v^o ' ^ 

S~aZ* *•*-» 


aF?\ 

'AdLf dt:J„_dS 0 


1 dL * 


where use has been made of the coefficients 
C( 1; 1) = 1 

C(2; 2) = 1 C(2;l, 1)=* 

C( 3; 3) = 1 <7(3; 2, 1)=# C(3; 1, 2) = * C( 3; 1, 1, 1)=*. 


5. ELIMINATION OF VARIABLES 

Since the solution of the system is known where F is reduced to F 0 , the problem is to eliminate 
variables which are not present in Fo. Suppose a canonical transformation is found in such a way that 
p of the n coordinates (p<n) have been eliminated from the Hamiltonian, that is 

F* - F*(l* + i, • • • A, LI ,Lt, •••,£:). (23) 

The equations of motion then yield 

L?-C*(const.) (fc- 1,2, •••,?). (24) 

If these constants are replaced in F*, then 

F*-F*(t . -,C„LU U " • ,Lt ) 

and the problem is reduced to one of n-p degrees of freedom, 
a. If p*n, the problem is completely solved, since 

Z*«C* (A— 1,2, •••,»») 

and 


(t *■ o>t(C i, • • • ,C,)H-f* (0) . (A «■ 1 ,2, • • • ,n) 
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b. If p = n — 1, the problem is integrable by quadrature. In fact, 

Lf = C k (& = 1, 2, 

t*^~A(C u C h £) 

Ui n 

«?■• c ” • • •> <v.,- L O- 

F*(Ci,"',C n -i ; (Z,LX) — C ~ const., 

L:=L;(c,c u --,c n - l ] o 

^=X,(C,C S ,-..,C^,;C) 

'C 


Since 

then 

and therefore 
and 


r & 


f) 


The coordinate (t becomes a known function of time as well as J Li. 


Therefore, the equations 


Ci, •••, C,. i; L*(<), C(0)(*=J. 2, n-1) 


can be integrated by quadratu..., 

The von Zeipel’s method consists in the elimination of some of the • > ordinates (angular variables) 
and the reduction of the problem to case (b) and possibly (a). T.n- adaptability of this method is 
based on a set of hypotheses which are listed below in Roman numeral::. 

(I) The new and old corresponding variables differ by a qii**-.' iit at least of the first order, i.e. 


t\~ti -0(e) 


(*-t, 2, 


Lt-L t ~ 0(«). 


This automatically fixes & to correspond to the identity transformation since for «-0, the above 
conditions give 

(fm ( t 


(t-1, 2 


LfmL t . 


Therefore 
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If expression (25) is substituted into Equations (19), (20), (21) and (22), then 


Fo(t, L*) = Ft((, L*) 


Lti — L/i Ct—li 

p =2 

pit r*\ . &F t . A/dSi\ dF<, 1 A /d»>i d$\ d*F 0 _ . r 

w - i)+ 4U A r *arr + 2 2^a?;j r T .eiw Fti,L j 

4 **±4 2fi4.lv d ' F * <9&\ 

+ 4 S dL i'f - 1 ~ dt * *,\ dL t)f . c dt; 5 * w *v 

and similarly for Equation (22). 

It is seen that S*(L*, ^) and F*(L*, l) are unknown functions. In order to perform a particular 
solution toward the elimination of certain angular variables in F* we impose conditions (which are usu- 
ally suitable in Celestial Mechanics) on the functions S* and . They are 
II) F*{L *, t) does not depend on /,(» = 1, 2,* • -,p<n) for any k> 0. 

Ill) S*(L*, t) only depends on the (,(i~ 1,2,* • *,n) through trigonometric functions, for • ny fc>0. 
This avoids “secular perturbations" in the momenta L )} or in other words differences 

T T tSp) 

ijJ L > bt, 

are periodic functions of the /*(&= 1, 2,* •*,«). 

The application of these conditions, together with the obvious fact that Po does not depend on angular 
varables £<(t = 1,2,* ♦ *,p<n) which are to be eliminated, yields the relations 

Po(W“, (n, (29' 

P „»P,* 


F | &Si\ &F o _ A / 6iS|\ bF* 

* bh t ^\bl*)^ ,**t 

<~pfi w 35 m 


P*+Pi,-P?+P& 


Fip+P% p + 


A/asA ifi P * . v (22l\ d lL 

**\dtjL,-LU L* 4 “ w d/, 


and so forth. The functions Pi, and Pi„ Fu and Pi„ P», and P»„ P*. and P*, are the portions of 
Pi, Pj, P* and Pi* which are respectively independent of and dependent on t>»e t t (\ «■ 1,2, • • \p),and where 

P v/^\ ^1.4.1 V * r * 

*" d^/, , *9lr2 h \ df t et, /, r 

it jm. 1 L*k l '*k 

p.. v ^ 1 . 1 V AMl *!l\ ni) 

■ 4 VatfA. , a/; 2 4. Vai-r St , at#; 1 

<-H-l ™ i i. J-P+l € k m ** 


life 
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In the usual problems of Celestial Mechanics F<> does not depend on any angular variable so that the 
Equations (30), (31), (32) and the corresponding equations for higher order are much simplified. Thus, 
the additional hypotheses will be considered. 

IV) F 0 and thus F* depend only on the momenta L* 

V) The angular variables ( t (i= 1,2, ••*,»») corresponding to momenta L t (i= 1,2,* • *,to) that are 
present in F„ have been eliminated to the k a order. 

The next problem is the possibility of elimination of angular variables whose conjugate momenta 
are not present in F 0 . At this stage the Hamiltonian of the system is 


where 


F* = F 0 *(L?,...,L:)+Ff(/« +1 ,..., /„ Lf,--., L*) + ...+Ff(/ m+1 ,..., 


L*=Cj=const(j=l, 2 ,♦••,*»), 


(33) 


and the old and new variables are related by 


r T * dSi , dSt , , 3Sk 




0 = 1, 2, ••*, n) 


• * dL * ' dL* 


dS/c 

‘ dL *' 


(34) 


Assume a new canonical transformation from the variables (C+ 1 , • * - i, • • • ,L*) to the 

variables (C+i,* ♦ • ,CV£Ti,* • *,!£*) and let 


S* = 3*(t*+ U l*u 

l>e its generating function. Then, since L* = C } = const 0=1, 2,* • •, to), 

F?(LTV ••,!£*) -f?*(LTV* *,!£•)- const Lt*=L{ = C k = const (fc= 1, 2,. • -,to) 
Ff(« +1 ,....«:C,,C s ...NC m .L*i,,--.Xr)=n*(C + .,---.«:L:ii L**). 


(35) 

(36) 

(37) 


The last equation implies that the elimination of further variables is possible if and only if Ff does not 
depend on them. For in this case 


Ft (t* +P + „ • • ;CCuC» * • .,C.,LWi, • • -,«* ) = FT* (« +F r „ • • • ; «•+ i,< • • • ,L ** ) 


and 








2 ' v m+P+l> 


•> <:■> • 


x < 3S > 

<-*+»>+i 1 


which defines <Sf. It is important to note that in such a case Si will be defined by an equation involv- 
ing 2nd order terms ; these terms are therefore necessary to obtain first order “perturbations." This fact 
is exactly what happens in Brouwer’s theory on artificial satellites (1959), where 

a) The elimination of g* is possible because F* is independent of this variable. 

b) The development for long period perturbations (those of argument g*) needs the evaluation 
of 2nd order terms. 



s#** iWUWW 
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This type of reasoning can be carried on up to any order in exactly the same way. It may then 
happen that the elimination of a certain angular variable by obtaining S* requires the evaluation of 
terms of the fcth order. 

However if F* depends on the angular variables to be eliminated the problem cannot be solved 
unless it happens that the remaining system has one degree of freedom. For example, this Is the case 
of the perturbations on the motion of an artificial satellite by the moon. 

6. SMALL DIVISORS 

The case of critical inclination for the theory of artificial satellites of an oblate planet for which Pi 
is the dominant zonal harmonic and J&& —J\, is a well known example of the problem of small divisors. 
Here, only a particular aspect of the question is dealt with. Consider the solution of Equation (30) 
in the usual case where F* does not depend on the t , . The characteristic associated system is 

dt i _ dti _ dtp _ dSi /pn-, 

dFjTdfo ~ dF 0 ~ F\ v ' K ) 

dL* dL* dL* 

Should one of the partials dF 0 /dL* happen to be zero, the general solution would certainly be discon- 
tinuous since a “small divisor” is present. However this divisor is not exactly zero because the quan- 
tity dFo/dL* is evaluated to first order only. 

In the case of critical inclination it is necessary to take 


(S = lSo - l - € 1 ^/Sl/2 - t'CiSl _ l'**^ l <Ss/2"l~ * * *. 

However, in doing so the separation of “long periodic” and “secular” perturbations is lost. The in- 
tegration leads, in most cases, to elliptic integrals (Hori, 1960). 

The question of small divisors usually arises whenever the problem presents cases of libration as 
particular solutions. 

Another case to be mentioned is the resonance for an artificial satellite whose period is commensur- 
able with the period of rotation of the Earth when tesseral harmonics are included. Again, expansion 
in powers of can be used to solve the problem (Morando, 1962). 

Finally it is important to note that singularities in the Equations (39) reflect singular points in the 
hypersurface defined by the Hamiltonian of the system in a phase-space of 2 (n-p) dimensions if p 
variables have already been eliminated. 


7 . SUMMARY 

The general differential equations of the von Zeipel’s method have been given to any order. It is 
hoped that this will avoid tedious Taylor expansions if one needs to go to order higher than the second. 

At the same time, the brief discussion on tho applicability and a few pathological cases of the 
method, may give some guidance toward the solu tion of new problems. 
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